import numpy as np
import matplotlib.pyplot as plt

# 设置辣椒素浓度为20 ppm
concentration = 20  # 单位：ppm

# 生成频率数据（聚焦1.3 - 1.8 THz区间）
x = np.linspace(1.3, 1.8, 1000)

# 预定义特征峰位置
peak1_pos = 1.42   # 主峰1
peak2_pos = 1.67   # 主峰2
shoulder_pos = 1.35  # 肩峰位置

# 基线函数（模拟油品的基础吸收谱）
baseline = 0.03 + 0.03 * x - 0.01 * x**2

# 辣椒素的两个主要特征峰函数
peak1 = 0.1 * (concentration / 100) * np.exp(-50 * (x - peak1_pos)**2)
peak2 = 0.1 * (concentration / 100) * np.exp(-50 * (x - peak2_pos)**2)

# 添加肩峰
shoulder_peak = 0.05 * (concentration / 100) * np.exp(-30 * (x - shoulder_pos)**2)

# 原始总吸收谱
y_original = baseline + peak1 + peak2 + shoulder_peak

# 计算二阶导数
y_second_derivative = np.gradient(np.gradient(y_original))

# 绘制图形
plt.figure(figsize=(10, 6))

# 绘制原始吸收谱
plt.plot(x, y_original, '-', color='blue', alpha=0.5, label='Original Spectrum')

# 绘制二阶导数光谱
plt.plot(x, y_second_derivative + 0.05, '-', color='red', alpha=0.8, label='Second Derivative Spectrum')

# 添加特征峰位置的竖虚线
plt.axvline(x=peak1_pos, color='k', linestyle='--', linewidth=1)
plt.axvline(x=peak2_pos, color='k', linestyle='--', linewidth=1)
plt.axvline(x=shoulder_pos, color='k', linestyle='--', linewidth=1)

# 添加肩峰分离的标注
plt.annotate('Separated Shoulder Peak', xy=(shoulder_pos, -0.01), xytext=(1.38, 0.01),
             arrowprops=dict(facecolor='black', shrink=0.05))

# 图表设置
plt.xlabel('Frequency (THz)')
plt.ylabel('Absorption coefficient α(ν) (cm⁻¹)')
plt.title('Second Derivative of Absorption Spectrum (20 ppm Capsaicin)')
plt.legend(loc='upper left', fontsize=8)

# 显示图表
plt.show()